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We investigate the nature of phase transitions induced by Dzyaloshinskii-Moriya (DM) interac- 
tions on a classical pyrochlore antiferromagnet. For both symmetry-allowed antisymmetric exchange 
interactions, the macroscopic degeneracy due to geometrical frustration is relieved and a long-range 
magnetic order appears in the ground state. We find an Ising-like phase transition in the case of di- 
rect DM interaction which selects a doubly degenerate all-in all-out magnetic order. In the presence 
of indirect DM coupling, the magnet undergoes an XY-like phase transition into a state with broken 
Z§ symmetry. We show that the critical behavior in both cases is modified due to the constrained 
spin fluctuations in the correlated Coulomb phase. We also demonstrate an interesting order from 
disorder phenomenon where the system switches between two distinct types of broken Z§ symmetry. 



I. INTRODUCTION 

Geometrical frustration refers to the inability of spins 
to satisfy conflicting interactions simultaneously due to 
the lattice connectivity. 1 Magnets with strong geometri- 
cal frustration exhibit a variety of unusual ground states, 
elementary excitations, and phase transitions. One of 
the most intensively studied frustrated systems is the 
Heisenberg antiferromagnet on the pyrochlore lattice. 
Monte Carlo simulations with interactions restricted to 
nearest-neighbor spins showed that the magnet remains 
disordered even at temperatures well below the energy 
scale of exchange constant. 2 Contrary to the uncorrelated 
paramagnetic state, motions of spins in this so called 
Coulomb phase are subject to a set of local constraints. 
In the continuum approximation, these local constraints 
translate to the Gauss's law for a fictitious magnetic 
field. 3 - 5 

Many interesting properties of the Coulomb phase can 
be traced to the macroscopic degeneracy in the clas- 
sical ground state of the pyrochlore magnet. Indeed, 
minimization of the nearest-neighbor exchange interac- 
tion requires that the total spin on every tetrahedron be 
zero, which still leaves an extensive number of uncon- 
strained degrees of freedom. 2 This degeneracy in turn 
makes the magnet susceptible to small perturbations such 
as anisotropies, dipolar interaction, etc. In general, the 
frustration is relived and some sort of long-range order is 
selected by the dominant perturbations. 

Phase transitions taking place in the highly con- 
strained Coulomb phase also exhibit interesting, and 
sometimes peculiar, behavior. A recurring theme in ex- 
periments and numerical simulations has been the ap- 
pearance of first-order phase transitions, despite that a 
Landau-type analysis of the symmetry-breaking phase 
would otherwise predict a continuous one. Discon- 
tinuous phase transitions have been observed in py- 
rochlore antiferromagnet perturbed by magnetoelas- 
tic coupling, 6 ' 7 long-range dipolar interaction, 8 ' 9 ex- 
change anisotropy, 10,11 and further-neighbor exchange 
interactions. 12 ' 13 It is known that strong fluctuations of 



the order parameter could induce a first-order transition, 
a scenario frequently occurs when the number of com- 
ponents of the order parameter is greater than four. 14,15 
This is indeed the case for most of the phase transitions 
mentioned above. 

The dipolar spin correlations stemming from the strong 
local constraints also affect the nature of phase tran- 
sitions from the Coulomb phase. 3-5 As is well known, 
the long-range dipolar interaction can change the uni- 
versality class of a magnetic phase transition, modify- 
ing its critical behavior to that of a model with short- 
range interactions in higher spatial dimensions. 16,1 In a 
recent study on the effects of the exchange anisotropy, it 
was demonstrated that Neel ordering from the Coulomb 
phase indeed displays a critical behavior similar to that 
of a uniaxial dipolar antiferromagnet. 10 

As most broken symmetries observed in pyrochlore 
lattice are rather complex, it is desirable to study 
phase transitions which can be described by simple or- 
der parameters. In this paper, we investigate criti- 
cal behavior of magnetic ordering due to asymmetric 
exchange coupling, or the Dzyaloshinskii-Moriya (DM) 
interaction. 2 ' 21 The relevant order parameters have ei- 
ther a Zi or Zq symmetry, depending on the sign of the 
DM interaction. Specifically, we study the classical spin 
Hamiltonian: 

H = Jj^ s i • s i + J2 D « ' ( Si x S ^ ' w 

(*i> (*i> 

where the two terms represent the isotropic exchange 
and the DM interactions, respectively, between nearest- 
neighbor spins. In addition to being an interesting spin 
model in itself, recent studies have revealed the impor- 
tant role of DM interaction in the formation of magnetic 
spirals 22 ' 23 and the magnon Hall effect 24 observed in the 
pyrochlore magnets. 

By applying the so called Moriya rules, 21 the possible 
forms of the coupling vectors D^- were determined by 
Elhajal et al. up to a multiplicative constant. 25 Depend- 
ing on the sign of this constant, the corresponding DM 
terms are called direct and indirect DM interactions. The 
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FIG. 1: Phase diagram of Heisenberg antiferromagnet 
with Dzayloshinskii-Moriya interaction on pyrochlore lattice 
(shown in the inset). The ground state of direct DM inter- 
action (D > 0) consists of all-in all-out tetrahedra. The or- 
dered phase of indirect DM interaction (D < 0) contains two 
distinct ordered states with broken symmetry: noncopla- 
nar and orthogonal structures at high and low T, respectively. 
The solid and dashed lines denote second-order transition and 
crossover lines, respectively. 



ground states for both forms of the DM interaction were 
also obtained in Ref. 25. In particular, it was pointed out 
that indirect DM interaction leaves a continuously degen- 
erate ground state at the mean-field level. Monte Carlo 
simulations, however, showed that a six-fold degenerate 
ground state with coplanar spins was selected by thermal 
fluctuations. 25,26 Although their preliminary simulations 
indicated a continuous magnetic ordering for both types 
of DM interaction, the authors did not touch upon the 
nature of the phase transitions. 

In this paper we performed detailed finite-size scaling 
to study the critical behavior of model (1). The phase 
diagram obtained from our extensive Monte Carlo sim- 
ulations is shown in Fig. 1. Consistent with previous 
results, 25 we found that the direct DM term is mini- 
mized by an all-in all-out spin structure, whereas the 
indirect DM interaction selects orthogonal magnetic or- 
ders as T — > 0. By carefully characterizing the symmetry- 
breaking phase in the indirect DM model, we find a new 
intermediate regime characterized by a magnetic order 
with noncoplanar spins. A crossover line separates the 
two different types of broken Zg symmetry. 

Due to the Zi symmetry of the doubly degenerate all-in 
all-out spin structure, magnetic ordering induced by di- 
rect DM interaction (D > 0) is expected to be Ising-like. 
The critical exponents obtained from finite-size scaling 



indicate that the transition is indeed close to 3D Ising 
universality class. We show that the constrained spin 
fluctuations in the Coulomb phase result in noticeable 
deviation especially in exponents /3 and 7. 

The case of indirect DM interaction (D < 0) is more 
complicated. Upon lowering the temperature, the mag- 
net undergoes an XY-like transition into a phase with a 
broken Zq symmetry. We introduce a doublet order pa- 
rameter £ = (CxiCy) such that the six coplanar ground 
states correspond to the six corners 9$ = nir/3 of a 
hexagon in the ^-domain; here n is an integer. Interest- 
ingly, below the critical temperature the magnet enters 
an ordered state characterized by 9^ = (n + 1/2) n/3; 
the corresponding magnetic order consists of noncopla- 
nar spins. As temperature is further lowered the system 
gradually evolves toward the coplanar magnetic struc- 
ture. The crossover between the two ordering patterns 
does not break further symmetries. 

The continuous degeneracy of the noncoplanar ground 
states observed in Ref. 25 manifests itself as an 0(2) sym- 
metry of the doublet order parameter £. Consequently, 
both types of broken Z§ symmetry arc entirely due to the 
cntropic effect, or so called order-from-disorder. By ex- 
plicitly computing the magnon contribution to the free 
energy, we demonstrate that the above crossover phe- 
nomenon is indeed dictated by thermal fluctuations. 

The remainder of the paper is organized as follows. In 
Sec. II, we introduce three antiferromagnetic order pa- 
rameters to characterize the spin configurations in the 
ordered phase. The critical behavior of direct and in- 
direct DM interactions is discussed in Sees. Ill and IV, 
respectively. A quantitative analysis of the crossover be- 
tween the two distinct types of broken Z§ symmetry is 
presented in Sec. V. Finally we conclude with a discus- 
sion of these results in Sec. VI. 



II. GROUND-STATE MAGNETIC ORDERS 

From the symmetry viewpoint, the DM interaction 
is allowed on the pyrochlore lattice since the nearest- 
neighbor bonds are not cnetrosymmetric. 21 In Moriya's 
original theory, the DM term stems from the rela- 
tivistic spin-orbit interaction. 21 Regardless of its mi- 
croscopic origin, the high symmetry of the pyrochlore 
lattice places strict constraints on the direction of the 
coupling vectors. 25 For a bond oriented along the [110] 
direction, the vector points along [110]: Dy = 
(±D, tD,0)/V2. The value of D y on any other bond 
can then be found through symmetry transformations. 

The expression of the DM interaction could be fur- 
ther simplified using the antiferromagnetic order param- 
eters Li = (So + Si — S2 — S 3 )/4S', 22 and so on (see 
Fig. 2 for labeling of the spins), which measure the stag- 
gered magnetizations of a tetrahedron. For complete- 
ness we also define the ferromagnetic order parameter: 
M = ^j = nSj/4S'. In terms of these variables, Hamil- 
tonian (1) can be recast into a sum over tetrahedra 
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FIG. 2: Magnetic states of a tetrahedron with zero total spin 
M = 0: (a) all-out, (b) orthogonal, and (c) noncoplanar struc- 
tures. The three antiferromagnetic vectors Li, L2, and L3 
measure the staggered magnetizations between the two red, 
green, and blue bonds, respectively. 

H = where the Hamiltonian of a single tetrahe- 

dron is 

Hr = 8JS 2 \M\ 2 -4V2DS 2 (a- L 2 x L 3 

+ b- L 3 x Li +c - Li x L 2 ). (2) 

Minimization of the exchange energy requires M = on 
every tetrahedron, which still leaves a macroscopic de- 
generacy. The remaining degrees of freedom encoded in 
the three staggered vectors are determined by minimiza- 
tion of the DM term. 

Since the vanishing of M makes the three staggered 
magnetizations orthogonal to each other and imposes 
a constraint on their lengths: ^2%=i l-M = I) 22 we 
parameterize the antiferromagnetic order parameters as 
= (pifii, where the three form an orthogonal triad 
and (pi > denote their lengths. Denning a 'handedness' 
for the triad: x = ' ^2 x n 3 = ±1, the DM term 
becomes 

-xDS 2 (hi ■ a 0203 + n 2 -hfofa + n 3 •60i0 2 ). (3) 

Note that inversion of the triad —> —hi also changes 
the chirality \ ~~ * ~Xi reflecting the time-reversal sym- 
metry of the DM interaction. In the following, we focus 
on the case of a right-handed triad with \ = +1. 

For D > the DM interaction is minimized first by 
aligning the triad to the principal crystal axes such that 
fix = a, and so on. Minimization of the resultant expres- 
sion subjecting to the constraint X)i=i §1 = 1 yields a 
ground state with 0i = 2 = 03 = l/v3- The staggered 
magnetizations are given by 

L 1 = ( X /v / 3)a, L 2 = ( x A/3)b, L 3 = (x/>/3)6, (4) 

with % = ±1 corresponding to the all-out [Fig. 2(a)] and 
all-in structures, respectively. 

The magnet remains frustrated in the presence of the 
indirect DM interaction, despite a significantly reduced 
degeneracy. The frustration comes from the fact that 
absolute minimum of the DM energy with D < and 
\ = +1 requires the unit vectors ii^ be antiparallel to the 
corresponding principal directions, e.g. hi = —a. How- 
ever, such a complete anti-alignment is impossible for two 



right-handed triads. Numerical minimization yields two 
inequivalent classes of continuously degenerate ground 
states. The first set consists of orthogonal spins lying in 
the plane perpendicular to one of the three cubic axes. 
The energy of a coplanar state is invariant with respect to 
uniform rotations of spins about the corresponding cubic 
axes. For spins perpendicular to the c axis, the magnetic 
state of the tetrahedron is described by 

Li = (cos a + sin b) /\/2, 
L 2 = (sin a — cos b) /y/2, 
L 3 - 0. (5) 

where the angle describes the uniform rotation about 
c axis. As T — > 0, thermal fluctuations select coplanar 
states with = or 7r. 26 The orthogonal structure with 
spins pointing along the [110] and [110] directions (0 = 0) 
is shown in Fig. 2(b). 

The second set of ground states consists of noncoplanar 
spins described by order parameters 

Li = cos(0 — 7r/4) cos 9 a, 
L 2 = sin(0 - 7r/4) cos 9 b, 
L 3 = — sin 9 c, (6) 

where 9 = arctan (v^sinf?) and 9 is a continuous param- 
eter. A tetrahedron with 8 — corresponds to orthogo- 
nal spins lying in the ab plane. Starting from 9 — 0, one 
can reach the other two inequivalent orthogonal states at 
9 = ±7r/4, where the spins are perpendicular to a and b 
axes, respectively. As will be discussed in Sec. V, order 
by disorder at finite temperatures prefers the noncopla- 
nar spins with, e.g. 9 — ir/2, shown in Fig. 2(c). 

As first noted in Ref. 27, the coarsed-grained antiferro- 
magnetic order parameters are just the three components 
of a flux or polarization field introduced in Refs. 4 and 
5 to describe the correlated Coulomb phase. This corre- 
spondence thus implies that the coarsed-graind staggered 
fields also obey a dipolar correlation at large distances. 
A continuum theory of the antisymmetric pyrochlore an- 
tifcrromagnet is outlined in Appendix A. 

Since the DM term does not depend explicitly on M, a 
remarkable feature of Hamiltonian (2) is that its ground 
state is independent of the antiferromagnetic exchange J. 
However, spin fluctuations in the J = and J — > 00 
(Coulomb phase) limits are quite different. In the latter 
case, minimization of the J term requires spin fluctua- 
tions satisfy J^ieiE ~ 0, which in turn result in a 
dipolar spin correlation at large distances. In the follow- 
ing, we compare the universality class of magnetic order- 
ing in the two limiting cases and show that the critical 
behavior is indeed modified by the constrained fluctua- 
tions in the Coulomb phase. 
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III. ISING-LIKE PHASE TRANSITIONS: 
DIRECT DM INTERACTION 

The direct DM interaction removes the degeneracy 
completely and selects a doubly degenerate all-in all-out 
ground state. The magnet undergoes a continuous phase 
transition at a temperature T c ~ O(D) as demonstrated 
by the specific heat curves shown in Fig. 3(a). The peak 
of the specific heat becomes sharper with increasing sys- 
tem size. However, finite-size scaling of the specific heat 
is rather difficult due to a nonzero regular component. 
Instead, by plotting the maximum of heat capacity vs 
the linear size of the system, we found a power-law de- 
pendence C max = Co + const x L a ' v , indicating a contin- 
uous phase transition. Here we have used an exponent 
ajv = 0.175 corresponding to the 3D Ising universality. 
To further characterize the phase transition, we turn to 
the finite-size scaling of the order parameter. 

In the all-in all-out structure, spins point along the 
local (111) directions with tetrahedra of opposite orien- 
tations in the all-in and all-out states, respectively. We 
introduce an Ising order parameter to characterize this 
doubly degenerate state: 



NS ^ ' 



(7) 



where N is the number of spins and points along the 
local [111] axes, i.e. from the center of type-I tetrahedron 
to the corresponding corner. The magnetic susceptibility 
and Binder's cumulant are defined accordingly 
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The temperature dependence of these variables are shown 
in Fig. 4 for various linear size L. As can be seen from 
Fig. 4(a), the order parameter \m\ starts to grow below 
a critical temperature T c rs 0.38J estimated from the 
crossing of the Binder's cumulant shown in Fig. 4(e). 




FIG. 3: Monte Carlo simulations for coupling constants 
D = O.UJ. (a) Specific heat C = ({E 2 } - {E} 2 ) /NT 2 as 
a function of temperature T for various system sizes, (b) The 
maximum of the specific heat C ma x as a function of L a '", 
where a/v = 0.175. The number of spins N = 16L 3 for a 
system with linear size L. 




FIG. 4: The temperature dependence of (a) order parameter 
\m\, (c) susceptibility \ m , and (e) Binder's cumulant B4 m 
obtained from Monte Carlo simulations for the specific case 
of D = O.UJ. The right panels (b), (d), and (f) show the 
corresponding finite-size scaling plots. The estimated critical 
exponents are j3 — 0.38, 7 = 1.121 and v = 0.63. 



The finite-size scaling plots shown in Figs. 4(b), (d), 
and (f) indicate that the data points indeed collapse on 
a universal curve. The critical exponents of the phase 
transition can be estimated from the corresponding scal- 
ing relations. For example, scaling of Binder's cumu- 
lant B 4m = B((T - T C )L X I V ) yields the critical expo- 
nent v = 0.63. Combined with the ratio ajv obtained 
from the specific heat, we find a = 0.11. Similarly, from 
the scaling of |m| and Xmi we obtain critical exponents 
P = 0.38 and 7 = 1.121. 

The discrete Z2 symmetry of order parameter m im- 
plies a phase transition in the 3D Ising universality class. 
We found that the critical exponents a and v indeed 
agree with the expected values for Ising transitions in 
three dimensions. The exponents /3 and 7, however, dif- 
fer from the corresponding values of the 3D Ising class 
[fi = 0.32 and 7 = 1.24). 28 As discussed in Sec. I, this dis- 
crepancy could be caused by the dipolar-like effective in- 
teraction between spins in the Coulomb phase. Nonethe- 
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FIG. 5: Monte Carlo simulations for coupling constants 
D = -0.7J. (a) Specific heat C = ({E 2 } - {E} 2 ) /NT 2 as 
a function of temperature T. (b) The maximum of the spe- 
cific heat C m ax as a function of L a/ ", where a/v = -0.0217. 



less, the obtained exponents roughly satisfy the identity 
a + 2(3 + 7 « 2. 

In Appendix B, we numerically study the magnetic or- 
dering of Hamiltonian (1) assuming J = 0. As discussed 
in the previous section, the same all-in all-out magnetic 
order is selected as the ground state. Finite-size scal- 
ing analysis with the order parameter m shows that the 
magnetic phase transition indeed belongs to the 3D Ising 
universality class in the J = limit. This contrasting 
result demonstrates that the discrepancy observed in ex- 
ponents /3 and 7 is attributed to the constraints imposed 
by a large J on the motion of spins in the Coulomb phase. 



IV. XY-LIKE PHASE TRANSITION: INDIRECT 
DM INTERACTION 

The indirect DM interaction only partially removes the 
magnetic frustration. At the single tetrahedron level, the 
degenerate ground-state manifold is parameterized by ei- 
ther (j> or 9 introduced in Eqs. (5) and (6), respectively. 
Nonetheless, the system undergoes a continuous phase 
transition at a temperature T c ~ £)(|D|). A clear peak 
can be seen in the temperature dependence of the spe- 
cific heat as shown in Fig. 5(a). To characterize the phase 
transition, we need an order parameter which transforms 
as an irreducible representation of the symmetry group 
and properly describes the ground-state manifold. 

To this end, we note that a general ground state of 
the indirect DM model can be constructed using the 
three orthogonal structures in which spins point along 
the (110) directions. We thus introduce unit vectors 
m§ )3 = (t1,±1,0)/V§ and mf )2 = ( T 1, T 1,0)/V2 for 
coplanar spins perpendicular to the c axis. Similar vec- 
tors m™ and can be defined for the other two coplanar 
states. For a given spin configuration {S^}, the projec- 
tion onto these orthogonal structures is given by 

p » = lvcX) Si ' ril i'' U> = a,b,c. (9) 
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FIG. 6: The domain of the order parameter £ = (Cx,Cy) i s 
bounded by a circle whose circumference corresponds to non- 
coplanar ground states described by Eq. (6). The coplanar 
ground states given in Eq. (5) lies on the three straight lines 
C, y = ±\/3Ce and £ y = 0. Also shown are the six orthogonal 
magnetic structures where spins point along (110) directions; 
The and ® symbols denote spin component coming out of 
and into the plane, respectively. They correspond to angles 
6q — rwr/3 where n = 0, • • • , 5. 

An order parameter q = max^p 2 was introduced in 
Ref. 26 to characterize the proximity of spin configura- 
tions to one of the six orthogonal structures shown in 
Fig. 6. Their Monte Carlo simulations demonstrated that 
q indeed approaches its maximum q ma _ x — 1 as T — > 0, 
indicating that the ground state consists of tetrahedra 
with orthogonal spins. 

As already mentioned in Sec. I, the ordered phase has 
two distinct types of broken Zq symmetry. To character- 
ize the different symmetry-breaking patterns, we intro- 
duce a two-component order parameter: 

Cx = ( Pa +Pb~2Pc)/V6, ( y = (Pa-Pb)/V2, (10) 

which transforms as the doublet representation of the 
group Oh- Note that the three projections are not in- 
dependent as their symmetric sum p a + Pb + Pc vanishes 
identically. The domain of possible values of the vec- 
tor £ = (£ x , <^y) is bounded by a circle with radius \/3/2 
[Fig. 6] . Its circumference is made of noncoplanar ground 
states which are parameterized by 6 [Eq. (6)]: 

\/3cos# „ 3sin# . . 

Cx = ~ =, {„ = -, =■ (11) 

V 2 + 4 sin 2 9 V 2 + 4 sin 2 9 

Noting that Cy/Cx — tan6^ = \/3tan(9, the six orthogo- 
nal structures shown in Fig. 6 correspond to 9^ = nir/3 
on the circumference, where n = 0, • • • 5. These six states 
are the magnetic ground state at T — >• 0. 

Another set of degenerate ground states is given by the 
three straight lines ( y — ±V3{, and Cy = 0. They cor- 
respond to the coplanar spins described in Eq. (5); the 
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FIG. 7: The temperature dependence of (a) order parameter 
C = | CI) ( c ) susceptibility XCi an d (e) Binder's cumulant B^q 
obtained from Monte Carlo simulations for the specific case 
of D = -0.7 J. The right panels (b), (d), and (f) show the 
corresponding finite-size scaling plots. The estimated critical 
exponents are /3 = 0.42, 7 = 1.181, and v = 0.672. 



Binder's cumulants at a temperature T c w 0.452 J indi- 
cates that the phase transition is second-order. 

The order parameter £ possesses a symmetry similar 
to that of XY model with a Z§ anisotropy term. It is 
known that this anisotropy perturbation is dangerously 
irrelevant in 3D and the ordering transition is in the 3D 
XY universality class. 29-32 To investigate the nature of 
the observed phase transition, we performed finite-size 
scaling on the relevant quantities. The critical exponent 
v which characterizes the scaling of correlation length can 
be estimated from the finite-size scaling plot of Binder's 
cumulant. As Fig. 7(f) shows, using v — 0.672 of the XY 
universality gives a rather good data collapsing. On the 
other hand, it is known that the XY model has a small 
negative exponent a = —0.0146. The nonzero regular 
component of the specific heat makes its finite-size scaling 
a rather difficult task. Fig. 5(b) shows the maximum 
of the specific heat as a function of scaled system size: 
C max = C + const x L a / V . Using a/v = -0.0217 from 
the XY universality class yields an agreeable result for 
large L. 

The finite-size scaling plot of the order parameter 
[Fig. 7(b)] shows that the data points indeed fall on a 
universal curve ( = L-PI V <$>({T - T C )L 1 / V ), from which 
we estimate the critical exponent j3 — 0.42. Simi- 
larly, scaling analysis of the susceptibility using relation 
Xc = L-'>/' / T{{T - T C )L X I") yields 7 = 1.181. These 
two exponents differ from the expected values of 3D XY 
universality class (/J = 0.3485 and 7 = 1.3177). 33 As al- 
ready discussed in the case of direct DM interaction, the 
discrepancy could be due to the constrained fluctuations 
of spins in the Coulomb phase. The long-range dipolar 
correlation of spins in this phase could modify the ex- 
pected XY critical behavior. Nonetheless, the deviations 
in /3 and 7 roughly compensate each other such that the 
critical exponents still satisfy the Rushbrooke equality 
a + 2/3 + 7-2. 28 



order parameter of an orthogonal structure parameter- 
ized by 4> is given by 

(x = ■\/3/2cos77 Al cos</>, ( y = \/3/2 sin 77^ cos <j>, (12) 

where f] a lh — ±27r/3 and rj c — 0, and <f> describes the 
global rotation of the coplanar spins about the corre- 
sponding cubic axes. 

Fig. 7(a) shows the order parameter as a function 
of temperature obtained from Monte Carlo simulations. 
The amplitude of the order parameter increases 

markedly below a temperature T c ps 0.45 J. In Figs. 7(c) 
and (e) , we show the temperature dependence of the sus- 
ceptibility and Binder's cumulant defined as 

xc = y «Kl a > - <ICI> a ) , ^c = i-^^. (is) 

Again, a clear peak which diverges with the system size 
can be seen in the susceptibility curves. The crossing of 



V. ORDER-BY-DISORDER AND CROSSOVER 
PHENOMENON 

A. Crossover between two types of Za order 

Although the order parameter £ is indeed nonzero in 
the ordered phase, the pattern of the broken symmetry 
is unclear from the above Monte Carlo simulations. Ex- 
pecting a state with broken Zq symmetry, we define an- 
other order parameter which is sensitive to the angular 
distribution of the vector £: 34 ' 35 

C 6 = (|C|cos(60 c )), (14) 

where 6^ = arctan(C l/ /Ca:)- The coplanar ground states 
with 9q = nir/3 is characterized by a positive ^6- I n P ar_ 
ticular, a maximum C6,max = \/3/2 indicates that the 
system is in one of the six orthogonal structures shown 
in Fig. 6. We refer to this symmetry-breaking pattern 
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as the type-I order. Another type of broken Z§ symme- 
try, referred to as the type-II order, corresponds to order 
parameter £ clustered around the six angular positions 
9^ = (n + 1/2) 7r/3. Tetrahedra in the type-II order con- 
sist of noncoplanar spins shown in Fig. 2(c). The Z$ 
order parameter is negative (q < for the type-II order. 

Fig. 8 shows the temperature dependence of the 
two order parameters and C = (I CI) obtained from 
Monte Carlo simulations. Depending on the angular dis- 
tribution of the doublet vector £ (shown in Fig. 9), we 
divide the ordered phase into three distinct regimes. In 
the XY regime just below T c , the distribution of the or- 
der parameter exhibits an emergent rotational symme- 
try [Fig. 9(b)]. This result is consistent with numerical 
studies on the closely related 3D six-clock model which 
confirmed the irrelevance of the anisotropy at T c . 

Surprisingly, instead of directly settling in the type-I 
order which is shown to be the ground state as T — > in 
Ref. 26, the system first enters a regime with type-II order 
as indicated by a negative (g. Upon further lowering the 
temperature, the system gradually switches to the type-I 
order. A crossover temperature can be estimated from 
the point where the order parameter ^ changes sign. 
At T — > the parameter approaches its maximum 
indicating a ground state with orthogonal spins. 

As discussed in Sec. IV, the continuous degeneracy 9 
of the noncoplanar ground state corresponds to a rota- 
tional symmetry of the doublet vector £. As a result, the 
observed broken Zq symmetry, either the type-I or type- 
II order, must be of entropic origin. Naively, one would 
expect quantum or thermal fluctuations select the six or- 
thogonal structures shown in Fig. 6 because they possess 
two zero modes corresponding to directions 9 and c6 in the 
ground-state manifold. By considering zero-point energy 
X)m £ ™/2 °f magnons, Canals et al. showed that quan- 
tum fluctuations indeed favor the six orthogonal struc- 



FIG. 9: Distribution of the instantaneous order parameter 
on the (CejCw) plane for (a) disordered phase, (b) rotation- 
ally symmetric regime, (c) type-II and (d) type-I phases with 
broken Zq symmetry. 

tures. However, at finite temperatures one must consider 
the magnon contribution to the system free energy. As 
we demonstrate below, the high-energy magnons actually 
play an important role in the high-temperature regime of 
the ordered phase. 

B. Order from disorder: Holestein-Primakoff 
expansion 

Here we compute the free energy due to the harmonic 
fluctuations around a given noncoplanar magnetic order. 
The free energy is a function of 9 which parameterizes 
the q = noncoplanar spins given in Eq. (6). We first 
introduce a local reference frame defined by three or- 
thonormal vectors: ef , ef, and n^, which are functions 
of 9. A small deviation from the noncoplanar structure 
can be expressed using the Holstein-Primakoff expansion: 

Si = s ( 1 - li) Ai + E < ^ + °^)- ( 15 ) 

a— x,y 

Here er^ = (erf, erf) whose components denote fluctua- 
tions along the two orthogonal local axes. Substituting 
the above expression into Eq. (1), we obtain a magnon 
Hamiltonian 

n = (j-V2D)J2\^\ 2 + zZ E H #°?°s> ( 16 ) 

i i^j a,b—x,y 

where 

fl# = i(Je?.e$+Dtf.e?x6$). (17) 



8 




FIG. 10: Dispersion of the lowest two energy bands near the 
r point. The gapless band contains a zero mode at k = 
corresponding to the direction 9 in the ground state manifold 
[c.f. Fig. 6]. The second band is gapped with A ~ D9 2 . 



Since the noncoplanar structure is translationally invari- 
ant (q = 0), the above Hamiltonian is diagonalized by 
Fourier transformation. In the vicinity of the coplanar 
state 0£ = 9 = 0, the lowest energy band has a gapless 
dispersion 



1 



£i,k« — (2J + 5V2Z>)|kj 



(4J + V2D)k 2 z , (18) 



where |kx| = kf. + fc 2 . The zero mode at k = cor- 
responds to the direction 9 in the ground-state manifold 
[see Fig. 6] . The dispersion of the second lowest band is 



1 



e 2 ,k « A + -(2J + V2D) |kj 



V2 " 



(19) 



where the energy gap is due to a finite DM coefficient 
and is proportional to the deviation 9 from the coplanar 
state: 



A w AV2D9 2 . 



(20) 



The vanishing of the gap at 9 — corresponds to the zero 
mode described by the <j> direction in Fig. 6, i.e. uniform 
rotation of spins about the cubic axes. At low temper- 
atures, spin fluctuations are dominated by the low en- 
ergy magnons. Due to the presence of two zero modes at 
9 = <fi = 0, the six orthogonal structures shown in Fig. 6 
are selected as the ground state as T — > 0, consistent 
with the Monte Carlo simulations. 

At higher temperatures, one needs to take into ac- 
count the contributions from the high energy bands. The 
magnon contribution to the free energy is given by 



^magnon^c) = £ Tlnsinh [e„, k (0 C )/T] . 



(21) 



n,k 



The angle 9^ = arctan(Cj//Cx) is related to the parameter 
9 through tan#^ = V3tan0. We numerically compute 




T=0.3J 



T=0.035J 



T=0.01J 



FIG. 11: Magnon free energy Eq. (21) as a function of angle 
6q = arctan(^ H /^). Note that this angle is related to the 
parameter 9 characterizing a general noncoplanar structure 
[Eq. (6)] through tan#<; = \/3tan#. The coupling constant 
used in the calculation is D = —0.2 J. Note that 9( = 7r/6 cor- 
responds to noncoplanar structure shown in Fig. 2(c), whereas 
9q — and n/3 correspond to orthogonal spins shown in 
Figs. 6 and 2(b). 



the free energy by summing over states within the Bril- 
louin zone. Fig. 11 shows the free energy as a function of 
0£, At higher temperatures, e.g. T — 0.3 J, the minimum 
of F magnoI1 occurs at 9q = (n + l/2)7r/3 corresponding 
to noncoplanar spins shown in Fig. 2(c). As the tem- 
perature decreases, local minimum starts to develop at 
angles 9q = nn/3 which correspond to the six orthogonal 
structures shown in Fig. 6 or Fig. 2(b). As temperature 
is further lowered, the orthogonal magnetic order takes 
over and becomes the absolute minimum of the free en- 
ergy. This calculation clearly shows that thermal fluc- 
tuations favor two distinct spin structures at different 
temperatures of the ordered phase, as indeed observed in 
Monte Carlo simulations. 



VI. CONCLUSION 

To summarize, we have studied the critical be- 
havior of a classical pyrochlore antiferromagnet with 
Dzyaloshinskii-Moriya interactions. Symmetry consid- 
erations allow for two distinct forms of antisymmetric 
coupling on the pyrochlore lattice which are termed di- 
rect and indirect DM interactions. 25,26 The macroscopic 
degeneracy of the nearest-neighbor Heisenberg spins is 
lifted for both types of DM interaction. Through exten- 
sive Monte Carlo simulations, we have demonstrated the 
continuous nature of magnetic ordering in both cases. 
The critical exponents of the phase transitions were ob- 
tained from finite-size scaling analysis. Finally, we have 
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shown that the ordered phase in the indirect DM model 
exhibits two distinct types of broken Z§ symmetry. 

The direct DM interaction completely relieves the geo- 
metrical frustration and selects a doubly degenerate all- 
in all-out magnetic order. We have observed an Ising- 
like phase transition with critical exponents a = 0.11, 
/3 = 0.38, 7 = 1.121, and v = 0.63. The deviation 
from the expected 3D Ising universality is attributed to 
the correlated spin fluctuations in the Coulomb phase. 
This is because minimization of the nearest-neighbor ex- 
change requires the total spin on every tetrahedra be 
zero, which in turn results in a dipolar-like spin correla- 
tion at large distances. As is well known, the long-range 
dipolar spin interaction could modify the nature of the 
magnetic phase transition. 

We have also simulated a pyrochlore magnet where 
the only coupling between spins is the nearest-neighbor 
DM interaction. Noting that this model selects exactly 
the same all-in all-out ground state, we have shown that 
the corresponding magnetic transition indeed belongs to 
the 3D Ising universality class when spin fluctuations are 
not subject to local constraints imposed by the nearest- 
neighbor exchange. 

To characterize the symmetry-breaking phase of the 
indirect DM model, we have introduced a doublet or- 
der parameter £. The domain of this doublet vector is 
bounded by a circle whose circumference is composed 
of degenerate noncoplanar ground states. On the other 
hand, coplanar ground states with spins perpendicular 
to one of the cubic axes lie on the lines ( y = ±V3( X 
and ( y = 0. Finite-size scaling analysis using the dou- 
blet order parameter yields a set of critical exponents: 
a = -0.0146, /3 = 0.42, 7 = 1.181, and v = 0.672. 
Similar to the case of direct DM interaction, we have 
observed noticeable deviations in (3 and 7 from the ex- 
pected 3D XY universality class. Again, the dipolar spin 
correlation in the Coulomb phase can cause the observed 
discrepancy. 

Our simulations have also uncovered an unusual or- 
der from disorder phenomenon in the case of indirect 
DM interaction. Thermal fluctuations below the transi- 
tion temperature lift the continuous degeneracy and se- 
lect a q = magnetic order with six-fold degeneracy. 
Interestingly, the magnet first develops a noncoplanar 
magnetic order before settling in the coplanar ground 
state. The crossover between these two types of Zq order 
corresponds to a rotation of the doublet vector £ from 
0£ = (n + l/2)7r/3 to 0q — mr/3 7 where n is an integer. 
By computing the magnon contribution to the free en- 
ergy, we have explicitly demonstrated that the observed 
crossover indeed originates from the entropic selection. 
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Appendix A: Continuum approximation 

Here we outline a continuum description of the py- 
rochlore antiferromagnet with DM interaction. Follow- 
ing Refs. 4 and 5, we define a magnetic or polarization 
field at the center of a tetrahedron as 



(Al) 



where ej denotes the local (111) direction at site i. Com- 
paring with the definition of the three staggered magne- 
tizations, we observe 



B a = 4S(L$,L%,L%) 



(A2) 



This relation was first noted in Rcf. 27. In the coarse- 
grained approximation, the constraint M = translates 
to V • B a (r) = for each component a. Expressed in 
terms of the staggered fields, we have 



d x L\ + d v L a 2 + d z L a 3 = 0. 



(A3) 



Noting that the states with small values of B a are entrop- 
ically favored, the probability distribution of the flux field 
has a Gaussian distribution, i.e. p cx e" H2a with 



(A4) 



The stiffness k is the single parameter of the theory which 
controls the amplitude of the correlations. The correla- 
tors of the staggered magnetizations at large distances 
are 



(L?(r)L$(0))«^ r ^ r6 3r<rj . 



(A5) 



While 7i2a is entirely of entropic origin, an additonal en- 
ergy term comes from the DM interaction 

n 2b = -V J d 3 r[a-L 2 (r) xL 3 (r) (A6) 

+b-L 3 (r) x Li(r)+c-Li(r) x L 2 (r)], 

where the coupling constant T> cx DS 2 . The competition 
between the two terms leads to a phase transition into an 
ordered phase with nonzero (Lj(r)). Details of the stag- 
gered fields are discussed in Sec. II for the two different 
forms of DM interaction. 

To ensure stability of the ordered phase and give penal- 
tics to short- wavelength fluctuations, we add gradient 
and quartic terms 

n 2c = J rf 3 r^[yl 1 |V-L l (r)| 2 + A 2 |VxL 4 (r)| 2 ],(A7) 



I thank E. Choi, N. Perkins, C. Fennie, and O. Tch- 
ernyshyov for collaboration on related works. 



H 4 



J d 3 r [uJ2 Mr)| 4 + vJ2 Mr)| 2 |Iy(r)| 2 ] . (A8) 

i ijij 
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The sum H2 + H4 thus constitute a continuum frame- 
work for describing the phase transition induced by DM 
interactions. Although the resultant expression resem- 
bles a conventional Landau-Ginzburg- Wilson functional 
expressed in terms of the antiferromagnetic order param- 
eters, it is important to note that these order parameter 
fields are not indepenent. Instead, they are subject to 
the constraint (A3). This condition can be satisfied by 
introducing a vector field A a = (Af, A% , A|) for each spin 
component a such that B a = V x A a . Regrouping the 
vector potentials A; = (Af,AV,Af), the coarse-grained 
staggered fields is given by 



Lj(r) = e ijk djA k (r) 



(A9) 



A detailed analysis of the resulting energy functional in 
terms of A^(r) will be presented in future publications. 



Appendix B: Ising transition in direct DM model 
with J = 

In this appendix, we examine the critical behavior 
of classical Heisenberg spins with direct DM interaction 
alone on the pyrochlore lattice: 



n = 



5> 



(Si 



(Bl) 



This model corresponds to the special case J = of the 
Hamiltonian (1) studied in the main text. Remarkably, 
the anisotropic direct DM interaction alone selects the 
same all-in all-out magnetic ground state [Fig. 2(a)] as 
in the J —¥ 00 limit. However, since spin fluctuations are 
not subject to constraints X^eKl ^» = imposed by a 
large J, this model provides us an opportunity to exam- 
ine how the universality class of the magnetic transition 
is modified by the long-range dipolar correlation of spins 
in the Coulomb phase. 

Fig. 12(a) shows the specific heat as a function of tem- 
perature. A clear peak which diverges with increasing 



(b) 






0' 




a' 








ps' 











1.2 1.3 1.4 1.5 



FIG. 12: (a) Specific heat C as a function of temperature T 
for various system sizes, (b) The maximum of the specific 
heat C m ax as a function of L a '", where a/v = 0.175. The 
number of spins N = 16L [i for a system with linear size L. 




FIG. 13: The temperature dependence of (a) order parameter 
\m\, (c) susceptibility Xm, and (e) Binder's cumulant B^m 
obtained from Monte Carlo simulations. The right panels 
(b), (d), and (f) show the corresponding finite-size scaling 
plots with 3D Ising critical exponents: /3 = 0.32, 7 = 1.24 
and v — 0.63. 



system sizes indicates a continuous phase transition at 
T c ~ O(D). Due to the existence of nonzero regular 
component in the specific heat, finite-size scaling of C 
is rather difficult. Fig. 12(b) shows the maximum of 
specific-heat as a function of L a ' v , where we have used 
the ratio a/v — 0.175 from 3D Ising universality class. 
The data collapsing on a straight line indicates a scal- 
ing relation C max = Cq + const x L a ' v . We use the 
same order parameter m defined in Eq. (7) to character- 
ize the doubly degenerate all-in all-out magnetic order. 
The temperature dependence of m and the correspond- 
ing susceptibility and Binder's cumulant are shown in 
Fig. 13. Using critical exponents of 3D Ising universal- 
ity: f3 = 0.32, 7 = 1.24, and v — 0.63, we obtained 
excellent data collapsing as shown in the corresponding 
finite-size scaling plots [Fig. 13]. The analysis shows that 
without the constraints imposed by nearest-neighbor ex- 
change J, the phase transition into the doubly degenerate 
all-in all-out state indeed belongs to 3D Ising class. 
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